Phase diagrams of classical spin fluids: the influence of an external magnetic field 

on the liquid-gas transition 
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The influence of an external magnetic field on the liquid-gas phase transition in Ising, XY, and 
Heisenberg spin fluid models is studied using a modified mean field theory and Gibbs ensemble 
Monte Carlo simulations. It is demonstrated that the theory is able to reproduce quantitatively 
all characteristic features of the field dependence of the critical temperature Tc{H) for all the three 
models. These features include a monotonic decrease of Tc with rising H in the case of the Ising 
fluid as well as a more complicated nonmonotonic behavior for the XY and Heisenberg models. The 
nonmonotonicity consists in a decrease of Tc with increasing H at weak external fields, an increase of 
Tc with rising H in the strong field regime, and the existence of a minimum in Tc{H) at intermediate 
values of H. Analytical expressions for Tc{H) in the large field limit are presented as well. The 
magnetic para-ferro phase transition is also considered in simulations and described within the mean 
field theory. 
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I. INTRODUCTION 

The investigation of continuum fluid models with cou- 
pled translational and spin degrees of freedom is of cur- 
rent theoretical interest [1-3]. The importance of such 
models lies in their property to display a rich variety 
of transitions between solid, liquid, and g well as 

magnetic ordered and disordered phases, which may oc- 
cur in real systems. The liquid-gas and magnetic para- 
ferro phase transitions in spin fluids were studied pre- 
viously by the mean field theory [3-7], the method of 
integral equations [8-13], and Monte Carlo (MC) simu- 
lation techniques [2,8,11,14-18]. The theoretical studies 
dealt mainly with spatially one-(d = 1) and thiee-{d = 3) 
dimensional Ising (n = 1) as well as three-dimensional 
Heisenberg (n = 3) fluids (here n denotes the spin dimen- 
sionality). In computer experiment, the magnetic tran- 
sition was investigated for the Heisenberg model [2,14] 
as well as for three- [16] and two-{d = 2) [17] dimen- 
sional Ising (n = 1) fluids using canonical MC simula- 
tions. The combined canonical and Gibbs ensemble MC 
(GEMC) simulations were performed for a Heisenberg 
system (d = 3, n — 3) to determine both the magnetic 
and liquid-gas transitions [8,11,15]. 

In most of the previous works, the liquid-gas coexis- 
tence was evaluated in the absence {H — 0) of an external 
magnetic field. Only a few papers [7,11,12] were devoted 
to theoretical study of the fluid behavior at iJ 7^ 0. It 
was found that for systems of hard spheres carrying Ising 
spins, an external magnetic field decreases the tempera- 
ture Tc of the gas- liquid critical point [7]. On the other 
hand, the presence of Heisenberg spins can lead to the 
inverse effect at strong enough fields [11,12]. As a result, 
a nonmonotonic behavior of Tc may arise in Heisenberg 



fluids due to a subtle interplay between the translational 
and spin degrees of freedom [12]. To our knowledge, no 
confirmation of the nonmonotonicity in Tc{H) has been 
given within computer experiment for continuum fluid 
models of spin systems. In the only work [11] done on 
GEMC simulations at 7? 7^ for a Heisenberg fluid, it has 
been concluded that the application of an external field 
increases the gas-liquid critical temperature Tc. But this 
conclusion has been made on the basis of results corre- 
sponding to just one finite (sufficiently large) value of H. 

It is worth mentioning also that, as far as we are aware, 
no theoretical calculations and computer simulations of 
Tc{H) have been performed for the planar {n — 2, d = 3) 
XY spin fluid model and no simulations on liquid-gas co- 
existence in the presence of an external field have been 
reported for the three-dimensional Ising fluid. Note that 
we are considering genuine fluid models in which spa- 
tial positions of spins are changed continuously (contrary 
to simplified so-called lattice gas schemes [19-21] where 
spins are allowed to occupy only positions belonging to 
sites of a chosen lattice). Moreover, all the works dealt 
with nonmagnetic repulsion interactions in the form of 
the simplest hard-sphere potential exclusively. In ad- 
dition, the magnetic interactions were truncated, as a 
rule, at some finite interatomic separation, without tak- 
ing into account long range corrections. The question of 
how these restrictions impact the behavior of Tc{H) has 
not been considered as well. 

In this paper we present a comprehensive study of the 
influence of an external magnetic fleld on the liquid-gas 
coexistence properties of fluid models with Ising, XY, and 
Heisenberg spin interactions. The corresponding phase 
diagrams are calculated in the whole region of varying H 
with the help of the GEMC simulation technique and a 
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modified version of tfie mean field theory. As is shown, 
a good agreement between the MC data and theoreti- 
cal results for Tc can be achieved for all the models and 
for any value H of the external field, including the limit 
i? — > 00. At i? = 0, it is demonstrated that the mean 
field theory predicts a tricritical point independent on the 
dimension n of the magnetic order parameter. The simu- 
lation results for the Ising model agree with this topology, 
whereas for the XY and Hcisenbcrg fluids a critical end 
point, beside a gas-liquid critical point in the ferromag- 
netic phase, is indicated. 



II. MEAN FIELD THEORY 
A. Models 

Let us consider three models of magnetic fluids for 
spatial dimension = 3 with spin interactions of Ising 
(n = 1), planar XY (« = 2), and Heisenberg (rt = 3) 
types. Within all these models, the total potential en- 
ergy of the A''-particle system can be cast in the form 



^ N N 



(1) 



Yi = {ri^,riy,ri^) denotes the spatial coordi- 
Si is the n-dimensional spin vector [i.e., (si^, 0,0), 



where 
nate, 

{si^,Siy,0), or {si^,Siy,Si^) for n = 1, 2, or 3, respec- 
tively] of unit length (|si| = 1), and = |ri — rj|. For 
convenience, the homogeneous external magnetic field 
H = (iJ, 0,0) is directed along axis X of the laboratory 
system of coordinates. The Yukawa function 



J{r) = — exp 
r 



a — r 



(2) 



is used to describe the internal magnetic interactions, 
where (T and e relate to the size of particles and coupling 
constant, respectively. The nonmagnetic interaction ip 
between particles can be modeled by the hard-sphere 
(HS) 



or soft-core (SC) 




(3) 



4e 



(7)"- (7)' 



e. 



0, 



< V2a, 



(4) 



repulsion potentials. Note that the attraction between 
particles is formed exclusively due to ferromagnetic in- 
teractions. This corresponds to a so-called "ideal" class 
of spin fluids, where the attractive part of nonmagnetic 
interactions is absent (see comments at the end of the 
paper). 



B. Equations of state 



Following the spirit of works [6,22], the Gibbs free en- 
ergy per particle corresponding to Hamiltonian (1) can 
be presented as 



(ln/(s)) 
/3 



JO 



+ / AaTj^ - H • 



(5) 



Here is the free energy of the reference system (in 
the absence of magnetic interactions and external fields), 
/(s) relates to the normalized single-particle hmction de- 
scribing the distribution of spins in orientational space, 
= k^T is the temperature with fee being the Boltz- 
mann's constant, ( ) denotes the statistical averaging, 
and 



m 



1 / ^ 

^ i=l 



s/(s)ds 



(6) 



defines the magnetization of the system. The contribu- 
tion to the free energy caused by spin interactions with 
the parameterized Yukawa function J a ir) can be written 
in the form 



j9a{r,Si,S2) 



dJgjr) 
da 



sr S2/(si)/(s2)drdsids2, 
(7) 



where p = N/V is the number density with V being the 
volume, and 5'a('', si, S2) introduces the pair distribution 

function of a system corresponding to a-switched on mag- 
netic interactions with Jq, = at a = and Jq. = J(r) 
at a = 1. 

Equation (5) formally leads to exact results but re- 
quires the knowledge of function gair, Si, S2) for each in- 
termediate states < a < 1. Since, in general, this 
function cannot be determined exactly, some approxima- 
tions are needed to calculate J^. Within the mean field 
(MF) approximation it is assumed that ga{r, si, S2) does 
not depend on a and accepts the form of the pair distri- 
bution function (r) of the reference system. The latter 
function can further be approximated by its values in the 
low density regime, g(p{r) « exp[—/3(p{r)]. Then in view 
of Eqs. (2) and (6), the integrations in (5) and (7) can 
be performed explicitly. This results in 



JF = + fceT y" ds/(s) In /(s) - ^ 



apm 



H-m, (8) 



where 



f° 

a = Att 

Jo 



g^{r)J{r)r^dv = 8^{T)'Kea^ 



(9) 



2 



is the magnetic interaction strength, and the multipher 
Jq°° exp[—p(p{r)]J{r)r^dr 



7(T) 



/~ J(r)r2dr 



(10) 



takes into account the softness of nonmagnetic repulsion 
potential (f. 

Considering the free energy (8) as a functional of /(s), 
it can be shown that the minimum of is achieved at 



where 



exp(/3h • s) 
J exp(/3h • s) ds 

h = H + apm 



(11) 
(12) 



can be treated as the effective field, consisting of the ex- 
ternal term H and averaged internal contribution apm. 
Eq. (11) defines, therefore, the cquifibrium single-particle 
distribution function in the MF approximation. Then, 
taking into account the fact that the vector m is parallel 
to H, the right-hand side of Eq. (6) can be integrated in 
quadratures. The result is 



'tanh(^). 



m = < 



H-\-apm 



T ( H+apm \ 

^0 \^ A;bT ) 



n = 1 , 

n = 2, (13) 



coth 



/ H+apn 



H+apm ' 



where Ii{x) = ^ e^'^°^^ cos{ltjj)d'ilj denotes the modi- 
fied Bessel function of the first kind and order I. Relation 
(13) represents the magnetic equation of state (MES) of 
the system. Note that the form of this equation depends 
on the number n of components of the magnetic order 
parameter m. 

The pressure equation of state (PES) can readily be ob- 
tained by partially differentiating Eq. (8) with respect to 
p, using the thermodynamic relation P = p'^{dT /dp)T,H 
with Eqs. (11) and (12). As a consequence, one finds that 
the total pressure is the sum of two terms, 



P 



Pip — —ap^rn^ . 



(14) 



namely, the pressure P^ = [p^dTip/ dp)T.H correspond- 
ing to the reference system and the part coming from the 
magnetization. For the HS reference system (3) we use 
the quasiexact Carnahan-Starling relation [23] 



P^{p, T) = pfcBT(l + 7? + r?' - rf){l - r?)" 



(15) 



with r? = -Kpa^ /& = 7rp*/6 being the packing fraction. In 

the case of a SC potential (4), the softness of Lp is taken 
into account by replacing in Eq. (15) the HS diameter 
a by its SC counterpart a^p. The latter quantity can be 



determined by requiring the second virial coefficients re- 
lated to the HS system with the particle's diameter cr^ 
and the SC system with the real potential (f to be equal 
between themselves. This leads to 



(16) 



where 



= ^J^ (l - exp[-/3<p(r)]) r^dr . (17) 

Then the SC pressure can be obtained using Eq. (15) 
with 7] = rj,p = irpa^/G = ^{T)np*/6. This is justified 
by the fact that the SC potential (4) is close enough to 
the HS function (3) (</'sc(^) increases rapidly to infinity 
with decreasing r in the range r < a, whereas it quickly 
tends to zero at r > a). 

The relations (9), (13) and (14) constitute the main re- 
sults of the MF theory. In the case of the HS potential (4) 
(when g^{r) = 1 for r > a and g^{r) — a,t r < a and, 
thus, 7 = ^ = 1, see Eqs. (10) and (17)), they coincide 
completely with those obtained earlier [5,6] (for n = 1 
and 3). Our expressions are more general, since they do 
not restrict us to the HS convention only, but also are di- 
rectly applicable for more realistic SC magnetic systems 
(including the case n = 2). 



C. Phase separations 

Analyzing the MES (13) at = 0, it can be shown 
that nontrivial (nonzero) solutions in the magnetic order- 
ing parameter m take place for temperatures lower than 
the Curie temperature T\ = ap/ {nk-o). In the dimension- 
less representation p* = po^ and T* = k^T/e, the mag- 
netic phase transition curve reads = 8TTp*^/{T^)/n. 
Since, in general, the function j{T) may depend on tem- 
perature in a characteristic way, the last equality repre- 
sents a nonlinear equation which should be solved with 
respect to at fixed p*. In the case of SC potential 
(4), the computations show that the deviations of 7(7") 
from unity do not exceed about 3% in a wide temper- 
ature range of 0.3 < T* < 6. For this reason we can 
put 7(T) = 1 without loss of precision (at least in the 
range mentioned above). Hence, the Curie temperature 
is found analytically, = Swp* /n. It linearly depends 
on the density and is inversely proportional to the num- 
ber of spin components. Near the Curie line at T < T^, 
the MES can be expanded in a series with respect to the 
deviation t={Tx- T)/Tx, yielding 



lim 



Tx 



Tx~T 



= c. 



(18) 



where c is a constant depending on the spin dimension- 
ality n (see Table I). 

In order to get the liquid-gas critical point T^., one has 
to look where the inverse compressibility goes to zero. 
The liquid-gas phase transition occurs on the Curie line. 
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so that in view of Eqs. (14), (15), and (18) one finds at 

T = Tx that 



/dP 



T,H 



d_ 
dr] 



1 



V- 



(1 - 77)3 



-f=0. (19) 



Thus solving Eq. (19) with respect to rj yields a tricritical 
point at pI = Qrit/n and = 8npl/n = ASijt/n. Within 
the HS version of the MF theory, the solutions % allow 
to be presented analytically: = 1 + — {432/ — 
72 - qy/^]/6 with q = 6(4i7)i/-^ - 24 - 96(2 /p)^/^ and 
p = 83 + 3^/993 forr^ = l,?7t = l + ^-il/^-l-q)^/^ 
with q = (pi/3 - 2 - 17/pi/3)/6 and p = 82 + 3\/T293 for 
n = 2, andr?t = 1 + [2^- (72/^- 24 - 4g)i/2]/6 with 
q = (p/2)i/3 - 26(2/p)i/3 - 2 and p = 245 + 9VT609 for 
n = 3. The values for and with three digits after 
the decimal point are collected in Table I and marked 



by the superscript HS. In the case of the SC version, we 
should replace (after taking the partial derivative) the 
packing fraction rj entering in Eq. (19) by its effective 
value 1]^ = ^(r)7rp*/6. It can be shown that for the SC 
potential defined by Eq. (4). the multiplier ^^/"^(T) de- 
creases monotonically from 1.05 to 0.94 with increasing 
the temperature in the interval T* G [0.3,6]. This be- 
havior has a simple physical meaning, namely, with in- 
creasing T the particle can approach to one another more 
closely due to the increase of their thermal velocities. 
As a consequence, the effective diameter = ^^^^{T)a 
will decrease. In such a situation (when ^ = ^(T*) with 
T* = 8TTp*/n), Eq. (19) transforms into a complicated 
nonlinear equation in p* and must be solved numerically 
(the integration (17) has been carried out numerically 
too). The results of these computations are shown in 
Table I as well and marked by the superscript SC. 



TABLE I. Mean field theory results in d = 3 



Model 


n 


c 


n 


„*HS 
Pt 




„*sc 
Pt 


^* SC 


Ising 


1 


3 




0.098 


2.462 


0.104 


2.621 


XY 


2 


2 


A'Kp* 


0.169 


2.119 


0.176 


2.207 


Heisenberg 


3 


5 
3 


Stt * 
TP 


0.224 


1.876 


0.229 


1.920 



In the presence of an external field (i.e., when H ^ 0) 
there is a liquid-gas phase transition curve ending in a 

critical point Tc. The critical temperature Tc and density 
Pc can be obtained numerically by solving the following 
system of two equations 



(20) 



where, in view of Eqs. (13) and (14), the pressure 
P{p,T,m) = P{p,T,m{p,T,H)) = P{p,T,H) should 
be considered as a function of p, T, and H. Therefore, 
the solutions Tc{H) and Pc{H) to system (20) must be 
found self-consistently with the solution m = m,{p, T, H) 
to nonlinear equation (13). The latter equation requires 
to be handled numerically as well. In the case n = 2, we 
have used a representation of the modified Bessel func- 
tions (appearing in Eq. (13)) in the form of the infinite 
series 7i(a;) = (a;/2)' X;^o(2=/2)^'=/(fc!('fe + 00 (restricted 
to a finite but large enough number of terms). 

The liquid-gas coexistence curve can be found at T < 
Tc by applying the Maxwell construction to pressure (14). 



Alternatively, we can introduce the chemical potential us- 
ing the relation p = J^+P/p and Eqs. (8), (11), and (14). 
Then one obtains 



p = p^- ksT In j exp(/3h • s) ds , (21) 



where 



exp(/3h • s) ds 



' 2cosh(ii±^, 

27r/o( 

sinh 



pm 



H+apm \ 



An- 



H+apm 



H+aprn 



n = 1 , 
n = 2, 

n = 3, 



(22) 

and /i^ is the chemical potential of the reference system 
which should take its Carnalian-Starling form 



p^ = k^T 



Inp-h 



?7(8 -9ri + 377^) 



(23) 
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to be sclf-consistcnt with Eq. (15). Note that the 
chemical potential and the pressure factor Z^{r],T) = 
Pip{p,T)/ (pkBT) are connected by the (exact) relation 
(3p^{v,T) = J^[Z^{v',T) - l]/v'dv' + Z^-1 + In(Ap) 
(the term feB^lnA, with A being the thermal de Broglie 
wavelength, has been excluded from the right-hand side 
of Eq. (23), since it depends only on T and is irrelevant 
for our consideration, see Eq. (24) below). The gas and 
liquid coexistence densities pg{T) and pl{T) are then de- 
termined applying the well-known mechanical and chem- 
ical equilibrium conditions 

P{PG,T) = P{puT) , /x(pG,T)=/z(pL,T). (24) 

In the regime of large magnetic fields, we can solve the 
MES (13) analytically taking into account the smallness 
of C = kB.T/{H -I- apm) «; 1. This gives 



l-2e '=bt +o{e ^) , n=l, 
l-5^ + ^(C')' " = 2, (25) 



H+ap 



where the terms of the second and higher orders for e~^/'> 
and ( have been omitted, and the inequality e~^/'» <C C 
has been used. Substituting (25) into the PES (14), tak- 
ing the derivative of P with respect to p, and solving the 
resulting equation {dP/dp)T,H = for T, we obtain the 
critical temperature as a function of H. The result is 



kBTcoo 



H ' 



n = 2. 



H 



n = 3, 



(26) 



where T^^o = lim/z^oo (-^ ) = ap^ool {k^W^^) and 
Pcoo = ^vcan^oo Pc{H) are the critical temperature and 
density in the infinite magnetic field limit, = 
d/dr][r]{l + r] + rf - rf)/{\ - r)f]L with -q^^ = 

I '/ — '/coo 

7rpcooC^/6, and it was assumed that the external field 
is much larger than both the kinetic energy and the in- 
ternal magnetic field, i.e., H ^ kBTcoo and H ^ apcoo- 
From Eq. (26) we conclude that with increasing the exter- 
nal field, the critical temperature approaches its limiting 
value from the top, when n = 1, and from the bottom, 
when n = 2 or 3. Note that the factor W'{ri) — 1 > is 
positive for all physical densities rj < 1. Moreover, for the 
Ising fluid model, the critical temperature Tc{H) tends 
to Tcoo exponentially with increasing H. This is not the 
case for XY and Heisenberg fluids, where Tc{H) reaches 
Tcoo slower, according to the inverse power law H~^. 



It is important to remark that in the saturation limit 
of infinite magnetic field H oo, all the spin fluid mod- 
els considered reduce to the same (nonmagnetic) fluid 
with the interparticle potential 4>{r) = ip{r) — J{r) con- 
sisting of the hard- or soft-core repulsion part as well 
as the Yukawa-like attraction. The reason is that then 
the spins align exactly along the field vector, so that the 
scalar product Si • Sj will be equal to 1 (see Eq. (1)) for 
any pairs of particles. The term H • will tends to a 

(infinite) constant and thus can be ignored in Eq. (1) (be- 
cause we are entitled to accept a new level for counting 
the energy of the system). The MF theory also leads to 
identical results for each n — 1, 2, and 3, when H ^ oo. 
Indeed, it follows from Eqs. (13) and (25) that the limit 
limH->oo TO = 1 is independent of n. The critical tem- 
perature and density at such a magnetic saturation, can 
be found as usually, using the general relation (20) with 
putting TO = 1 in Eq. (14). This yields p^^ w 0.249 and 
T*^^ « 2.264 as weU as w 0.262 and T*^'^ w 2.380 
for the HS and SC versions, respectively. Note also 
that in the limit H oo, expressions (21) and (22) 
for the chemical potential can be reduced (by extracting 
an infinite constant depending only on H) to the form 
p = Pip — ap that corresponds to a nonmagnetic system 
with the potential (j). 

Simulations show (see the next section) that for the 
Yukawa-fluid (YF) potential 4'{r) = (p{r) — J{r) with the 
SC repulsion, the critical liquid-gas temperature is equal 
to « 2.680. It is somewhat higher (within 12%) 
than the temperature T*^ k. 2.380 obtained within the 
SC version of the MF theory. It c;an be assumed that 
a signiflcant part of the above temperature discrepancy 
may come from the approximate form used for the equa- 
tion of state (sec Eq. (15), where rj = S^np* /6) of the SC 
reference system. Therefore, it becomes quite natural 
to introduce an adjustable SC (ASC) version of the MF 
theory, where the second virial parameter ^ (T) of the ref- 
erence system (Eq. (17)) is replaced by its rescaled analog 
|(T) = b^{T). The constant b can then be determined by 
requiring the critical temperature Typib) coincides with 
the exact result. This leads to 6 w 0.902 and corresponds 
to a slight decrease (by a factor of b^^^ « 0.966) of the ef- 
fective diameter when calculating the pressure according 
to Eq. (15) with T] = f(r)7rp76. 



III. COMPUTER SIMULATION AND THEORY 
CALCULATIONS 

A. Simulation procedures 

Two kinds of MC simulations have been carried out 
to investigate the critical behavior. First, in order to 
determine the liquid-gas coexistence properties we have 
applied the Gibbs ensemble MC (GEMC) approach [24] 
within an advanced biasing scheme [8] for handling spin 
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degrees of freedom. The GEMC simulations were per- 
formed for two system sizes, namely, with N = 500 and 
1000 particles (or only with N = 500 or 1000 for some 
temperatures and nonzero values of H to save computer 
time). In each case, the particles were distributed over 
two boxes with volumes which fluctuated under the con- 
straint of fixed total volume. The GEMC configurations 
were generated in cycles, where one cycle consisted of ei- 
ther (i) N trial displacements and spin reorientations of 
particles chosen at random; or (ii) one attempted volume 
rearrangement of the boxes; or (iii) A^ex attempts to ex- 
change the particles between the simulation boxes. The 
type of each cycle was selected at random and with equal 
probability among the above three possible steps. 

The bias has been used during the spin reorientation 
and exchange (insertion) steps, so that a new attempted 
direction for vector Sj was generated with a probabil- 
ity that favors orientations parallel to the local magnetic 
field. This bias was taken into account when consider- 
ing the acceptance probabilities corresponding to steps 
(i) and (iii) as well as when calculating the chemical po- 
tential. The analytical expressions for such probabili- 
ties have been obtained by modifying the hard-sphere 
acceptance rules [8] to the case when the soft- (instead of 
hard-) core nonmagnetic repulsion potential is present. 
The acceptance ratios for the particle moves and vol- 
ume changes were adjusted to lie in a range of 30 - 60%. 
The value for iVox was chosen to yield a success rate of 
particle transfers of 0.1 - 3%, depending on the temper- 
ature. The chemical potential was calculated during the 
exchange (insertion) step using the generalized Widom's 
method [25], while the pressure was obtained employing 
the virial theorem. 

Secondly, the Binder crossing technique [2,26] has been 
utilized to study the magnetic phase transitions (in the 
absence of an external magnetic field). Here, the usual 
canonical MC simulations have been performed for sys- 
tem sizes of iVi = 250 and = 500 particles at several 
particle number densities p and temperatures T. At each 
value of p, the Binder parameter B = 1 — M4/(3M|) 
was plotted as a function of T for the two system sizes 
(this parameter represents a fourth-order cumulant with 
Ml = (s') and s = |X]jSi|). The critical tempera- 
ture T\{p) of the para-ferro phase transition was then 
obtained from the position of the intersection point of 
curves Bn-^{T) and Bn^(T), i.e., from the condition 
BnAT\) = Bn^{Tx). As in the case of GEMC simu- 
lations, the orientational biasing technique [8] was used 
to improve the convergence of the canonical calculations. 

Within both the GEMC and canonical simulations, 
the Yukawa function was truncated at half of the box 
edge. In order to reduce the finite-size effects, a long 
range correction was taken into accoimt (within periodic 
boundary conditions) by including an additional term to 
the potential energy of the system. This term was de- 
rived [8] by integrating the Yukawa function beyond the 
cutoff radius and assuming that the binary distribution 
function is equal to unity in this region. For sufficiently 



large system sizes (as in our case) such an assumption 
should lead to virtually exact results. No additional fi- 
nite size scaling corrections were applied in determining 
the critical temperatures, assuming that they are small 
on the scale of phenomena considered. The number of 
MC cycles used to achieve an equilibrium state varied 
from 2 X 10^ to 5 X 10^ depending on the system and ex- 
ternal field value. Note that in the GEMC simulations, 
the pressure and chemical potential should be the same 
(within statistical noise) at equilibrium in the two (gas 
and liquid) boxes. After achieving the equilibrium, the 
investigated quantities were measured by averaging over 
2 X 10^ to 2.5 X 10^ cycles (depending on the model, en- 
semble, thermodynamic point, and field strength). The 
statistical uncertainties have been estimated using the 
block averages method [25] . 

The critical temperature Tc and density pc have been 
evaluated at ^ by fitting the discrete set of GEMC 
data near the criticality to the dependence p± = Pc + 
Ci(l - T/Tc) ± C2(l - T/Tc)f^ (the constants Ci and 
C2 should minimize the deviations between p± and sim- 
ulation values). Here the law of rectilinear diameters 
(pl + Pg)/2 = Pc + Ci{\ — T/Tc), the power law behav- 
ior pi^ — pQ = 26*2(1 — T/TcY , and the critical exponent 
P = 0.32 have been assumed as in Ref. [27] (where results 
for a Lennard- Jones fluid were analyzed). Note that for 
H the magnetic para-ferro phase transition disap- 
pears and the liquid-gas coexistence curves behave like 
those of nonmagnetic fluids with the effective attraction 
potential Jctf{r) = — (si • S2)HJ(r). At H = 0, where 
both the liquid-gas and magnetic para-ferro transitions 
exist simultaneously, the application of the above fitting 
should be performed with care because the rectilinear and 
power laws will not work properly at temperatures close 
enough to the temperature Tx of the magnetic transition. 



B. Results and discussion 

1. Ising fluid 

Examples of the liquid-gas coexistence obtained in the 
GEMC simulations for the soft-core Ising fluid at various 
values, H* = H/e ~ 0, 0.1, 0.5, 1, 5, and 00, of the ex- 
ternal field are shown in Fig. 1. As can be seen clearly, 
the critical temperature Tc goes down monotonously with 
rising the field strength and rapidly tends to its minimal 
value in the infinite field limit. For instance, already at 
H* = 5, the gas and liquid binodal branches are prac- 
tically indistinguishable from those corresponding to the 
case H 00. In addition, with increasing H the shape 
of the coexistence curves becomes wider near a critical 
point. At = 0, all the simulated points belonging to 
the gas phase lie on a curve which is very close in form to 
a straight line left to the magnetic para-ferro transition 
Curie line. The intersection of these lines defines a criti- 
cal end point which coincides with the liquid-gas critical 
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point expected for a tricritical point, i.e., Tec = = Tt- 
In other words, the Ising fluid exhibits a tricritical behav- 
ior. The shape of the Ising binodals has been verified by 
the multiple-histogram reweighting (MHR) method [28] 
and agreement with the GEMC results was observed. 



3.0 
2.9 

2.8 
2.7 
2.6 
2.5 
2.4 
2.3 



para- 




_ ferro — m&gnetic fluid, regions 




\h= 








H=QA \ 




k 


H=0.b 'A 




\ / 




gas 




H= CO W 

\\\ ferromagnetic 


o 




liquid 


" 1 


1 


metastable 


10 IT 




coexistence TJt 


■ s 1 

' £ i 

- s « 

- fi 

/ill 
"'1 


line 


region VA 


Curie 





0.0 0.2 0.4 0.6 0.8 
p* 

FIG. 1. The liquid-gas coexistence curves obtained from 
GEMC simulation data (circles) for a soft-core Ising fluid at 
various values of the external magnetic field H. The para- 
ferro phase transition (at -ff = 0) is plotted by the dashed 
(Curie) line passing via open squares (obtained in canonical 
MC simulations). The critical points are shown as crosses, 
whereas the tricritical point is presented by the star. The 
values of H are given in the dimensionless form = H/e. 



In making the above conclusion on the phase diagram 
topology at = we cannot be absolutely sure, since 
the GEMC simulations (as well as the MHR technique) 
do not provide us with precise enough data for temper- 
atures which are very close to the critical region. Near 
this region, the particle fluctuations become too large, so 
that the two GEMC boxes (which consist of finite num- 
bers of particles) can switch their identity many times 
during the simulations. This prevents one from obtain- 
ing good mean density values of the system in gas and 
liquid phases at T ~ Tc (although the density distri- 
butions over the boxes can still indicate the existence 
of two phases). On the other hand, having a discrete 
set of GEMC coexistence points lying relatively far from 
the criticality region, it is impossible to apply the fitting 
procedure (see at the end of the preceding subsection) 
in the vicinity of the expected tricritical point, where 
T ^ T\{pc) = Tc, because then the liquid and mag- 
netic transitions are coupled. The question concerning 
the topology of critical points at = goes beyond the 



scope of the present paper and requires additional in- 
vestigations by more sophisticated simulation techniques 
(such as finite-size scaling [29], for example). 

The primary goal of this work is to study the liquid-gas 
phase transition in the presence [H ^ 0) of an external 
field with focusing on the calculations of Tc as a function 
of H . The results of these calculations obtained in the 
case of the Ising fiuid within HS, SC, and ASC versions 
of the MF theory are shown in Fig. 2 together with the 
simulation data. Here, a more extended set of external 
field values has been used. It can be seen that the HSMF 
and SCMF schemes are able to predict qualitatively the 
monotonic decrease of the critical temperature with 
rising H. The relative deviations between the simulation 
data and SCMF predictions for Tc are of order 10 15%, 
i.e., they are not so large in view of assumptions made 
within the MF approach. When the ASC version is used, 
the theoretical and simulation data appear to be practi- 
cally indistinguishable and the disagreements do not ex- 
ceed the GEMC statistical noise. At weak external fields 
{H* < 0.5), the function Tc(if) behaves as T^{0) - cH^'^ 
(see Fig. 2, where the latter dependence is plotted at 
c = 0.258 for the ASC case). This confirms the MF pre- 
diction \imH^o[Tc{H)-T^{Q)] oc -iJ^/s^ obtained previ- 
ously in Refs. [7,30]. In the strong field regime {H* > 2), 
the critical temperature decreases exponentially, accord- 
ing to the analytical formula derived in sect. II. C (see the 
first line of Eq. (26)). 
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FIG. 2. The critical temperature Tc as a function of the 

external magnetic field H obtained for the soft-core Ising fluid 
from GEMC simulations (circles) in comparison with the re- 
sults of the HS, SC, and ASC versions of the MF theory (soUd 
curves). The critical level corresponding to the infinite value 
of the external field is plotted by the horizontal dashed line. 
Vertical bars indicate statistical uncertainties. 
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The Ising gas-liquid coexistence curves evaluated 
within the ASC version of the MF theory as well as 
the MF para-ferro transition line are presented in Fig. 3. 
The deviations in the theoretical and simulation binodals 
(please compare Figs. 1 and 3) are larger than the dis- 
crepancy in the case of function Tc{H) (sec Fig. 2) and 
they cannot be reduced to zero even within the ASCMF. 
For instance, the theory somewhat overestimates the gas 
phase densities and underestimates the densities in the 
liquid phase. As a result, the shape of the binodals ap- 
pears to be narrower with respect to that of the simula- 
tion coexistence curves. A more acciirate theory should 
be applied to describe quantitatively the liquid-gas phase 
diagrams in the whole density, temperature, and mag- 
netic field ranges. 



and Tc{H) presents, in fact, a rather complicated tech- 
nical problem. This is especially true in the case of XY 
system (n = 2) and SC/ASC versions of the MFT, where 
the integration in Eq. (17) as well as the computation of 
Bessel's functions are required additionally. Similar diffi- 
culties arise in MC simulations of XY fluids when apply- 
ing the orientational biasing technique [8]. The reason is 
that then the trial orientation vectors should be gener- 
ated with distributions, which cannot be presented (con- 
trary to the cases n = 1 and n = 3) by simple algebraic 
expressions, and the use of time-consuming Bessel-like 
functions is required. Nevertheless, developing an effi- 
cient algorithm has allowed us to overcome the technical 
difficulties and calculate the coexistence curves for all the 
models, including the XY. 
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FIG. 3. The gas-liquid coexistence of the Ising fluid ob- 
tained within the ASC version of the MF theory. The MF 
para-fcrro transition is plotted by the dashed line. For other 
notations see the caption to Fig. 1. 



2. XY-spin fluid model 

As was pointed out in the introduction, until now no 
computer experiment and theoretical investigations on 
the liquid-gas coexistence have been performed for the 
planar XY spin fluid model. In this respect it should 
be mentioned that in order to obtain the coexistence 
curves within the MF theory, one has to solve numeri- 
cally Eq. (24). At the same time, the HS, SC, and ASC 
functions Tc{H) are obtained by finding numerical so- 
lutions to Eq. (20) at given values of H. Taking into 
account that the magnetic equation of state (13) must 
be handled numerically too, the calculation of /9l,g(^) 
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FIG. 4. The gas-liquid (curves) and para-ferro (dashed line) 
phase coexistences obtained from GE (at different values of 
an external field H) and canonical (at H = 0) MC simula- 
tions, respectively, for the XY-spin fluid model. The critical 
end point is shown by the star. Other notations are the same 
as in Fig. 1. 



The XY gas-liquid binodals evaluated in the GEMC 
simulations are presented in Fig. 4 for the set H* = 0, 
2, 5, 10, 40, and oo of external field values. As can be 
clearly seen, here the topology of phase diagrams differs 
in several aspects from that of the Ising fluid. First, at 
H = the para-ferro magnetic line (which is included in 
Fig. 4 as well) intersects the gas branch of the binodal 
at a critical end point. Tee, which does not coincide with 
the critical point Tc, i.e., |Tc — Tce| ^ 0. So that, contrary 
to the Ising fluid, where the gas can only be in a para- 
magnetic state and the liquid be only ferromagnetic, the 
XY system exhibits a richer pattern. Here, the gas can 
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be either in paramagnetic (when T < Tec) or ferromag- 
netic (when Tec < T < Tc) states, and thus the hquid-gas 
transition can take place with keeping the ferromagnetic 
ordering. Secondly, with increasing H, the critical tem- 
perature starts to decrease rapidly reaching its minimal 
value at H* ^ 2 and further begins to increase much 
slower tending to the infinite- field limit Tcoo- The shape 
of the binodals near Tc{H) also demonstrates a nonmono- 
tonic behavior, namely, it first becomes wider with rising 
H and then begins to narrow, but the most sharp shape 
remains at H = 0. 
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FIG. 5. The critical temperature as a function of the ex- 
ternal field obtained for the XY-spin fluid model from GEMC 
simulations (circles) in comparison with the results of the HS, 
SC, and ASC versions of the MF theory (solid curves). Other 
notations are the same as in Fig. 2. 



The nonmonotonic dependence of Tc on the strength 
H of the external field observed from the simulations for 
the XY fluid is shown in more detail in Fig. 5. The corre- 
sponding results of the HS, SC, and ASC MF approaches 
are included there too. Again we can see that the HS and 
SC versions lead to qualitatively correct results and the 
deviations from the simulation data do not exceed about 
15%. In addition, all the approaches predict a minimum 
of Tc{H) at nearly the same external field H* k, 2. More- 
over, the ASC approach provides us with virtually exact 
results for Tc{H) at any value of H. For instance, the 
theoretical discrepancy is less than the level of GEMC 
uncertainties. However, at if ^ one has to be careful 
in application of the MF iJ^/^-dependence (see subsec- 
tion in.B.l) to interpretation of the function Tc{H) build 
on simulation data. Here, a crossover to an analytic be- 
havior with limij^o dTc/ dH = is expected because of 



the lack of a tricritical point (|Tc — Ted 7^ 0)- For large 
fields (H* > 50), the critical temperature tends to the 
infinite-field value Tcoo as Tc{H) — Tcoo oc —1/H, which 
is consistent with MF prediction (26) at n = 2. The sat- 
uration regime, when Tc{H) is almost equal to Tcoo, can 
be achieved here by very strong fields of order H* > 500 
(for the Ising fluid the saturation level is much lower, 
H* > 5, see Fig. 2). 
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FIG. 6. The gas-liquid and para-ferro coexistences of the 
XY-spin fluid observed within the ASC version of the MF 
theory. Notations are similar to those of Fig. 4. 



The theoretical liquid-gas coexistence curves and para- 
ferro transition line of the XY fluid are shown in Fig. 6. 

Comparing these results with the corresponding GEMC 
data (see Fig. 4) we see that, as in the case of the Ising 
system, the ASCMF predictions of pc and pl are not as 
perfect as those of the field dependency of Tc. Again, at 
H ^ the shape of the binodals becomes narrower be- 
cause of the overestimation and underestimation of gas 
and liquid densities, respectively. For H = 0, the AS- 
CMF approach does not predict the existence of critical 
and critical end points found in the simulations, but leads 
instead to a tricritical behavior. 



3. Heisenherg fluid 

The liquid-gas coexistence curves and the Curie line 
obtained for the Heisenberg fluid within the simulations 
and ASCMF theory are shown in Figs. 7 and 8, respec- 
tively. As can be seen, the Heisenberg system exhibits a 
topology of phase diagrams similar to the XY fluid. In 
particular, according to the simulation results iov H = Q, 
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the Curie line ends at a critical end point on the gas 
side of the binodal. so that Tec < 7c, where Tc should 
be referred to the liquid-gas critical point located in the 
ferromagnetic phase (see Fig. 7). It is worth mentioning 
that some evidence of the lack of a tricritical point in the 
Heisenberg system has been provided by GEMC simula- 
tions earlier [8]. Later on, it was stated [15] that owing to 
finite size effects it is very difficult to come definitely to 
one of the two possible scenarios: whether the Curie line 
ends at a critical end point (as suggested by the simula- 
tions, see Fig. 7), or at a tricritical point (as suggested by 
the MF theory, see Fig. 8). Relatively recently, the exis- 
tence of a critical end point for the Heisenberg fluid has 
been observed within the integral equation theory [12]. 
However, the difference jTce — Tc| reported there was very 
small and amounted to about 0.1%, which can be compa- 
rable with numerical uncertainties. On the other hand, 
this difference is much smaller than the discrepancy in- 
troduced by the mean spherical approximation, used in 
Ref. [12] for the closure relation to the Ornstein-Zernike 
equation. 



respectively. Again, the theory is able to qualitatively 
describe the liquid-gas coexistence properties (see Fig. 8), 
whereas the function Tc (H) can be calculated within the 
ASCMF approach quantitatively. The HS and SC ver- 
sions of the MF theory also reproduce well the critical 
temperature, but they underestimate Tc{H) to within an 
order of 10%. The minimum of Tc(iJ) at iJ* « 2 can 
be predicted by either version of the MFT. It should be 
pointed out also that the presence of SC repulsion leads 
to an increase of the gas-liquid critical temperature with 
respect to that of the HS potential (the SCMF curve lies 
above the HSMF at all values of i?, see Figs. 9 and 10). 
The same conclusion is valid for the Ising and XY models 
(see Figs. 2 and 5). The asymptotic behavior of Tc(7J) at 
iJ* » 1 looks like Tc{H) - Tcoo (x-l/H. It is identical 
in form to that of the XY system, but differs in the value 
of the coefficient of the proportionality to 1/H. Because 
this coefficient is twice as large at n = 3 (see Eq. (26)), 
the saturation regime will begin here at stronger fields 
{H* > 1000). 
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FIG. 7. The gas-liquid (curves) and para-ferro (deished 
line) phase coexistence obtained from GE (at different values 
of an external field H) and canonical (at H = 0) MC simu- 
lations, respectively, for the soft-core Heisenberg fluid. Other 
notations are the same as in Fig. 4. 



FIG. 8. The gas-liquid and paxa-ferro coexistences of the 
Heisenberg fluid obtained within the ASC version of the MF 
theory. Notations are the same as in Fig. 6. 



In the presence of an external field, the change in shape 
of the liquid-gas binodal and the change of the critical 
temperature are nonmonotonic for the Heisenberg fluid 
like for the XY system. This can be seen from Fig. 7 for 
the set H* = 0, 2, 5, 10, 15, 50, and oo, as well as from 
Figs. 9 and 10, where the values of Tc{H) are presented 
for a more complete set of H and in more detail at small 



Note also that in the limit H 0, the MF func- 
tion Tc(F) behaves like ~ -if^/s (^jth dT^dH ~ 
—H~^/^ — oo) independently of n (see section IV). 
However, as in the case of XY fluid, the function Tc{H) 
corresponding to simulation data for the Heisenberg 
model (see Figs. 9 and 10) should exhibit an analytic 
bahavior at ^ with limn^o dTc/dH = 0, because 
of the existence of a critical end point (|Tc — Tce\ ^ 0). 
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FIG. 9. The critical temperature Tc as a function of the ex- 
ternal field H obtained for the Hciscnbcrg fluid from GEMC 
simulations (circles) in comparison with the results of the HS, 
SC, and ASC versions of the MF theory (solid curves). Other 
notations are the same as in Fig. 5. 
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FIG. 10. The same as in Fig. 9, but the behavior of Tc(i/) 
at small values of H is shown in more detail. 



tation of all the functions using only one graph. As can 
be seen, the simulation points agree well with the AS- 
CMF predictions, although the deviations are somewhat 
larger than in the case of the field dependency Tc{H) 
of the critical temperature. The main part of these de- 
viations should be associated with MC statistical noise. 
For the absolute values of p*^ = limij^oo Pc' related to 
the ASCMF theory and GEMC simulations, we have ob- 
tained 0.295 and 0.307, respectively. The nonmonotonic- 
ity of pc{H) for n = 3 and 2 is closely connected with 
the corresponding nonmonotonous behavior of Tc{H). 
With switching on the external field and its slight in- 
creasing, the critical density begins to increase rapidly 
for all n = 1, 2, and 3. Further, the Heisenberg and XY 
dependencies Pc{H) exhibit a maximum at H* « 2, i.e., 
approximately at the same point, where the functions 
Tc{H) have a minimum (see Figs. 5 and 9). For n — 3 
and 2, after reaching the maximum, the functions Pc{H) 
decrease tending to the same limiting value pcoo- On the 
other hand, for n = I the function Pc{H) continues to 
increase monotonically to Pcoo 

on the whole interval of 
varying H (at the same time, the function Tc{H) mono- 
tonically decreases with rising H, see Fig. 2). 




FIG. 11. The normalized critical density Pc/pcoo as a func- 
tion of the applied external magnetic field H, evaluated within 
the ASCMF theory (curves) for the Ising (n = 1), XY (n = 2), 
and Heisenberg (n = 3) spin fluid models. The values ob- 
tained from GEMC simulations are presented by triangles, 
squares, and circles, for n = 1, 2, and 3, respectively. 



For completeness of our consideration we present in de- 
tail in Fig. 11 the normalized dependencies Pc{H) / Pcoo 
of the critical density on the external field for the Heisen- 
berg as well as XY, and Ising fluids, obtained within the 
ASCMF theory in comparison with the GEMC results. 
The normalization allows to make a quite visible presen- 



Finally, the simulation and theory results on the para- 
ferro coexistence are collected for n = 3, 2, and 1 in 
Fig. 12. It repeats to some extent the Curie lines already 
shown in Figs. 1, 3, 4, 6-8, but represents them in a con- 
siderably wider temperature interval. Despite the fact 
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that only a set of three MC points (for each n) is avail- 
able in a restricted density region, it can be stated that 
the dependence of Tx on p is almost linear and excellently 
coincides with the MF straight line T* = ^p* (the MC 
uncertainties are of order of the size of the symbols in 
Fig. 12). 
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FIG. 12. The paxarferro coexistence lines obtained within 
the MF theory for the Ising (n = 1), XY (n = 2), and Heisen- 
berg (n — 3) spin fluid models in comparison with canonical 
MC simulation data (circles). 



IV. CONCLUDING REMARKS 

In the present study we have obtained a complete set 
of phase diagrams for a class of Ising, XY, and Heisen- 
berg "ideal" spin fluids. The phase diagrams have been 
calculated applying a new version of the MF theory as 
well as the Gibbs- and canonical-ensemble MC simula- 
tion techniques. The new version takes into account the 
softness of nonmagnetic interparticle repulsion and cor- 
rects the equation of state of the reference system. This 
has allowed us to describe quantitatively the dependence 
of the critical temperature of the liquid-gas transition 
on the strength H of an external magnetic field for all 
the models considered over the whole region of varying 
H. Let us summarize and discuss the main results found: 

(i) It has been established that a common feature in- 
herent in the function Tc{H) is its decrease with turn- 
ing on and increasing H. In the case of the Ising model 
(n = 1), such a behavior remains valid for any further in- 
crease of H , including the infinite field regime [H oo). 
For the XY and Heisenberg models {n = 2 and 3), the 
decrease of with rising H gradually transforms into 



the inverse dependence, i.e., into an increase of Tc at in- 
termediate and strong fields. In this respect the discrete 
Ising model (where the spins can accept only two values, 
-hi or —1) exhibits a very specific feature which is not 
observed for the XY and Heisenberg models with contin- 
uous spin distributions. 

(ii) From a physical point of view, the field effects just 
mentioned are caused by the existence of two compet- 
ing mechanisms. In order to make the explanation more 
clear, one uses some expressions of the MF theory (which 
describes already the main features of the field depen- 
dency of Tc and pc). Then the spin fluid can be treated as 
a simple nonmagnetic system with an effective attraction 
potential Jes{r) = —(si -82)// J(r), where (si-S2)iy = 
within the MF approach. This corresponds to the pres- 
sure P = P^—acfip^/2 with Ues = am? being the effective 
attraction strength (see Eq. (14)). Taking into account 
the dependence of m on p (and H), one obtains that 
UbTc = paeff(l + o,px)/W{p)^ where x = dm/dH is the 
magnetic susceptibility of the system. This result follows 
from Eq. (15), the condition {dP/dp)T, = 0, and the MF 
relation dm/dp = amx (see Eq. (13)). Differentiating Tc 
with respect to H yields 



dH 



amp 



2x(l + apx) + amp 



dx_ 
dH 



(27) 



where the first and second terms on the right-hand side 
should be associated with the contributions of two dif- 
ferent mechanisms. The first mechanism is due to the 
fact that the external field favors to align the spins along 
H and thus increases m, i.e., x > for arbitrary val- 
ues of H. For genuine nonmagnetic systems, where Oeff 
is independent of p, one obtains the well-known result 
fceTc = pacs/W . This obviously means that the critical 
temperature is higher for fluids with stronger attractions 
between particles (e.g., the liquid-gas transition disap- 
pears completely for systems with only repulsive interac- 
tions, when fleff = 0). In our case, the effective attrac- 
tions aeff = am^ will rise with increasing the magnetic 
field because of rising m, and therefore the critical tem- 
perature will increase also (note that the term 2x(l-|-apx) 
remains always positive, see Eq. (27)). 

The second mechanism is more subtle and caused by 
the dependence of x on H . At strong enough fields 
{H 3> fceT and H ^ ap), when the magnetization m is 
very close to its saturation value, we obtain from Eq. (25) 
that X = dm/dH ~ exp[-2iJ/(fcBT)] forn = 1 as well 
as X ~ 1/^^^ for n = 2 and 3. Then dx/dH\n=i ~ 
-cxp[-2H/{kBT)] < and dx/dH\n=2,3 ~ -1/if^ < 0. 
Thus for n = 2 and 3, the influence of the second term 
in the right-hand side of Eq. (27) on the value dTc/dH 
can be neglected (this term is higher order of smallncss 
in H with respect to the first contribution). In other 
words, the first mechanism dominates and this indeed 
is observed for the XY and Heisenberg models, where 
Tc increases with rising H in the large magnetic field 
regime. This is, however, not the case for the Ising 
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fluid, where the second term (being negative) appears 
to be the same order in H and even greater in magni- 
tude than the first contribution, owing to the exponential 
field dependency. Then the critical temperature will de- 
crease with increasing H. In particular, for — > oo, 
combining the two terms gives limH^ca kBdTc/dH = 
-8{W^^ - l)exp[-2H/{kBTcoo)] < 0. 

At weak magnetic fields ~ 0, the second mechanism 
will prevail for all the models. Following arguments pre- 
sented in Ref. [7] at n = 1, it can be shown within the MF 
approximation that dT^/dH ~ -C(n)i7~^/^ < for any 
n = 1, 2, and 3, where C(n) > is the n-depcndent co- 
efficient of the proportionality. In order to better under- 
stand why the second mechanism suppresses the liquid- 
gas separation, one rewrites the term ampdx/dH in the 
equivalent form pdx/dp — apy^ . Then adding it to the 
first term 2x(l -|- apx) results in 2^(1 + \cipx) + pdx/dp. 
Thus, the negative contribution of the second term into 
the derivative dTc/dH directly follows from the fact that 
the susceptibility of the gaseous phase is larger than that 
of the coexistent dense liquid phase, i.e., dx/dp < 0. This 
means that the magnetization of the gaseous phase grows 
stronger, making the two phases more indistinguishable 
from one another. This effect is very strong in the limit 
H ^ 0, where the gas branch of the binodal is very close 
to the Curie line (along which x ^ oo). At the same 
time, the liquid branch quickly deviates from this para- 
ferro transition line, leading to especially large negative 
values of dx/dp (with limr^Tc dx/dp — oo a,t H ^ 0). 

(iii) We see, therefore, that for the Ising fluid, the sec- 
ond mechanism dominates (owing to the discrete nature 
of spin reorientations in this case) at all values of H. This 
explains the monotonic decrease of with rising H for 
n = 1. With increasing the number n of spin components 
to 2 and 3, the ability of the external field to decrease 
the critical temperature becomes smaller, and the first 
mechanism begins to prevail at larger values of H. This 
leads to the nonmonotonic behavior oiTc{H) for the XY 
(n = 2) and Hciscnberg (n = 3) models. The behavior of 
pc{H) can also be explained by appealing to the interplay 
between the same two competing processes. 

(iv) It has been shown that the MF theory does not 
predict at = the existence of a critical point and 
critical end point, found in the simulations for n = 2 
and 3, but leads instead to a tricritical behavior, ob- 
served in the simulations for n = 1. In this context it, 
first, should be pointed out that we dealt with a specific 
class of spin fluids, where the attraction between parti- 
cles is due to ferromagnetic interactions. According to 
the general classification, this corresponds to an infinite 
value of the ratio R = J g^{r)J{r)dr/ J g^{r)ifiattr{r)dr 
of the integrated strengths of the magnetic interactions 
and an attractive part of the nonmagnetic interactions, 
R = oo (because (fatti = 0). For finite values of R, 
different types of the phase diagram topology can be ob- 
served [5,7] within the MF theory (at i? = we come 
to the usual nonmagnetic fluids). Secondly, for R = oo 
and H = 0, we have concluded that the question on the 



topology of phase diagrams is very delicate, and requires 
more accurate theories and more sophisticated simulation 
techniques. These and related topics will be addressed in 
further investigations. 
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